import numpy as np
from scipy import signal
import matplotlib.pyplot as plt

# 1. 生成示例数据
t = np.linspace(0, 10, 1000)  # 时间轴
mean_velocity = 5.0          # 平均速度（m/s）
fluctuation = 1.0 * np.sin(2 * np.pi * 2 * t) + 0.5 * np.random.randn(len(t))  # 脉动分量
velocity = mean_velocity + fluctuation

# 2. 计算平均速度
avg_velocity = np.mean(velocity)

# 3. 提取脉动速度（减去平均速度）
fluctuating_velocity = velocity - avg_velocity

# 4. 去趋势处理（移除线性趋势）
detrended_fluctuation = signal.detrend(fluctuating_velocity)

# 5. 计算脉动速度的均方根（RMS）
rms_fluctuation = np.sqrt(np.mean(detrended_fluctuation**2))

# 6. 计算紊流度
turbulence_intensity = rms_fluctuation / avg_velocity

# 输出结果
print(f"平均速度: {avg_velocity:.3f} m/s")
print(f"脉动速度 RMS: {rms_fluctuation:.3f} m/s")
print(f"紊流度: {turbulence_intensity:.3f}")

# 可视化
plt.figure(figsize=(12, 8))

plt.subplot(3, 1, 1)
plt.plot(t, velocity, label='总速度')
plt.axhline(avg_velocity, color='r', linestyle='--', label='平均速度')
plt.ylabel('速度 (m/s)')
plt.legend()

plt.subplot(3, 1, 2)
plt.plot(t, detrended_fluctuation, label='去趋势脉动速度', color='orange')
plt.ylabel('速度 (m/s)')
plt.legend()

plt.subplot(3, 1, 3)
plt.plot(t, detrended_fluctuation**2, label='脉动速度平方', color='green')
plt.xlabel('时间 (s)')
plt.ylabel('平方值')
plt.legend()

plt.tight_layout()
plt.show()